Force networks and elasticity in granular silos 
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We have made experimental observations of the force networks within a two-dimensional granular 
silo similar to the classical system of Janssen. Models like that of Janssen predict that pressure within 
a silo saturates with depth as the result of vertical forces being redirected to the walls of the silo 
where they can then be carried by friction. By averaging ensembles of experimentally-obtained force 
networks in different ways, we compare the observed behavior with various predictions for granular 
silos. We identify several differences between the mean behavior in our system and that predicted 
by Janssen-like models: We find that the redirection parameter describing how the force network 
transfers vertical forces to the walls varies with depth. We find that changes in the preparation 
of the material can cause the pressure within the silo to either saturate or to continue building 
with depth. Most strikingly, we observe a non-linear response to overloads applied to the top of 
the material in the silo. For larger overloads we observe the previously reported "giant overshoot" 
effect where overload pressure decays only after an initial increase [G. Ovarlez et al., Phys. Rev. E 
67, 060302(R) (2003)]. For smaller overloads we find that additional pressure propagates to great 
depth. This effect depends on the particle stiffness, as given for instance by the Young's modulus, E, 
of the material from which the particles are made. Important measures include E, the unscreened 
hydrostatic pressure, and the applied load. These experiments suggest that when the load and 
the particle weight are comparable, particle elasticity acts to stabilize the force network, allowing 
non-linear network effects to be seen in the mean behavior. 

PACS numbers: 45.70.-n,45.70.Cc,83.80.Fg,45.05.+x 



I. INTRODUCTION 

The physics of granular matter is of great interest both 
because of the enormous breadth of physical systems in 
the granular regime as well as the lack of a characteri- 
zation of the granular state in fundamental terms. Al- 
though there has been some success understanding ener- 
getic granular gases in terms of kinematic models, dense 
granular materials continue to pose difhcult problems. 
We consider dry granular matter at low humidity where 
there are no attractive forces between the constituent 
particles and the particles do not significantly interact 
with the interstitial fluid (air). Contact forces dominate 
dense granular materials in this regime, with each parti- 
cle having several persistent contacts. 

Coulomb's law of friction is an inequality relating the 
magnitude of the friction force, Ff to the normal force, 
Fn, depending upon the mobilization of the friction. If 
/is is the static friction coefficient, then the mobilization 
is given by Ff /{Fnfis)- For a non-moving contact the 
mobilization ranges from —1 to 1 and reflects both the 
direction in which the force of friction acts to resist mo- 
tion and how close the forces acting on the contact are 
to causing motion in the opposite direction. Because the 
force of friction is dependent upon how the contact arose, 
the granular state is determined both by the positions of 
the particles and the history of how the particles came 



* Current address: National Center for Computational Toxicology, 
US EPA, Research Triangle Park, NC 27711 
t Electronic address: , wambaugh@phy. duke . edu | 



to be in these positions [ij . Determining a statistical de- 
scription of these contacts is at the core of understanding 
dense granular matter. 

Grain-scale statistical descriptions are further compli- 
cated by the presence of force networks |3, Isl]. An ex- 
traordinary feature of granular contact networks is that 
the distribution of forces on the contacts can be strongly 
inhomogeneous Q- Certain sequences of contacts carry 
forces with magnitude many times the mean for long dis- 
tances along "force chains" . As shown in Fig. [U force 
networks can develop to support the weight of the parti- 
cles as well as any additional applied load. 

Despite the difficulties inherent to understanding gran- 
ular materials, phenomenological engineering models ex- 
ist that describe the mean behavior well enough to allow 
for the design of granular storage facilities and handling 
processes. These models can provide useful tools for in- 
vestigating granular materials. In particular, Janssen's 
model of granular materials in vertical silos is interest- 
ing because it has been thoroughly tested, is realizable 
in the laboratory, and ultimately is viewed to be a work- 
ing qualitative and even quantitative description of the 
mean behavior of granular matter Q . 



II. THE JANSSEN MODEL 

Janssen treats the material in a silo as a continuum 
and considers horizontal slices, for which the difference 
of stress, the weight of the material in the slice, and 
the force of friction at the wall must be balanced ver- 
tically. Here we consider a quasi-two-dimensional silo, 
such as the one we use in the experiments reported be- 
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FIG. 2: Janssen's model predicts that stress exponentially 
approaches a saturation value with depth in a silo. Overloads 
decay exponentially and if an overload equal to the the satu- 
ration stress is placed on the silo the stress is constant. Here 
we plot predictions in a system that saturates at 0.5 in arbi- 
trary units for no load, a saturation overload {— 0.5), and an 
above saturation overload. 



FIG. 1: A force chain network in a silo filled with photoelastic 
disks under a 56 g load 



Wxy\ < fJ'wl'^xxl for wall friction /i^j. Janssen further as- 
sumes that the vertical and horizontal normal stresses 
are linearly related to each other to write axy = k^^ayy 
where k represents both the mobilization of interparti- 
cle contact friction (with friction coefhcient 6) and the 
mobilization of friction at the walls. By considering the 
extremes of maximum upward and downward mobiliza- 
tion, it has been shown that < k < 0. 

Using Janssen's relation for the shear stress in Eq. [1] 
allows us to rewrite it as a differential equation for ayy: 



da. 



yy 



dy 



(2) 



If Janssen's parameter k is taken as a constant, then this 
equation permits exponential solutions of the form: 



'yy 



oe 



(3) 



where $o is any load applied to the top of the system 
and 



low. Janssen assumes that the stress does not vary hor- 
izontally, allowing one to write for the case of a two- 
dimensional silo: 



A = 



L {(Jyy (y + A) - ayy [v)) = pffLA 



2(72; j^A 



(1) 



where x and y are the horizontal and vertical directions 
respectively, A is an infinitesimal change in y, a is the 
stress tensor, p is the density of material (mass per unit 
area), g is the acceleration of gravity and L is the width 
of the silo. The second term on the right-hand-side cor- 
responds to the two frictional wall contacts — one at 
each side of the silo. Coulomb's law of friction gives the 
shear stress at the wall in terms of the horizontal stress, 



2/c^^ 



(4) 



sets a length-scale for the evolution of stress with depth. 
As indicated by Fig. [21 the Janssen solution predicts that 
any overload will decay exponentially, and that the stress 
within the material will saturate with depth. This screen- 
ing effect results from friction at the walls acting to carry 
the weight of the material and any overload. Although 
we have presented the two-dimensional version here, the 
Janssen model is qualitatively similar in three dimen- 
sions. 

Janssen's prediction of stress saturating with depth 
has been validated qualitatively in numerous three- 
dimensional experiments that examine the apparent 
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weight of a confined granular pile. Though these experi- 
ments in effect only measured the forces at the bottom of 
the pile, by repeatedly varying the amount of material — 
and thus the height of the pile — a mean des crip tion of 
stress with depth was obtained 0, H, H, [13, [ill, lU ■ Large 
scale simulations of three-dimensional systems have also 
observed stress saturation [Tsj . In quasi- two-dimensional 
experimental systems that are similar to the analysis 
above, the saturation length scale has been found to scale 
with system size as predicted, A oc L [T3 |. 

While qualitatively accurate, Janssen's predictions 
have been shown to underestimate the actual value of 
saturated stress in a silo [l^ . Careful analysis of stress 
profiles in experiments and simulations have shown that 
for small depths the profile is more linear than exponen- 
tial [13, nil- These deviations from the Janssen model 
are not surprising, since key assumptions of the analysis 
can be shown to be incorrect. 

For the Janssen model above we assumed that the redi- 
rection parameter k was constant. For a cohesionless 
granular material, however, the ratio of the horizontal to 
vertical stress is constant only if the horizontal and ver- 
tical directions are the frame in which the stress tensor 
is diagonal [J. This is inconsistent, since Janssen as- 
sumes that there is a non-zero a^y at the walls. If axx 
and ayy are not the principal stresses, then we expect 
k ~ o'xx/o'yy to Vary horizontally and with depth. 

Several factors contribute to variation in the redirec- 
tion parameter with depth. In particular, it is difficult 
to achieve a uniform mobilization of friction at the walls. 
Experiments and simulations have shown that agreement 
with Janssen predictions can be obtained if the granular 
material is slowly lowered relative to the walls to make 
uniform the direction and extent of mobilization for all 
wall contacts [3, [l^, [l3] ■ The global packing fraction of 
the material has also been shown to alter the redirec- 
tion parameter and it is reasonable to expect that local 
variations iiipacking fraction would produce variations 
with depth [8|. Finally, experiments studying the stress 
within a silo of deformable beans, instead of the more 
rigid glass or metal balls typically used, have found that 
the redirection parameter varies with depth due to the 
elasticity of the particles [Tsj . 

More elaborate analysis of the Janssen model using 
stochastic differential equations to account for random- 
ness in the mobilization of friction has shown that the 
Janssen model o nly correctly predicts the stress on the 
axis of the silo [la]- In general, the vertical stress at 
a given depth has been found to vary horizontally, fur- 
ther complicating the force balance analysis [3, [3 ■ In 
three dimensions, attempts have been made to account 
for radial stress balance by using vertical force balance 
on concentric rings about the axis with the general result 
that force saturation due to screening by the walls is still 
expected, though the approach to the saturation value 
may not be exponential [20|. 

Although screening-induced saturation has been vali- 
dated, the behavior of overloads can be strikingly differ- 



ent from predictions. Instead of a monotonic, exponen- 
tial decay to the saturated stress supported by the walls, 
a non-monotonic, "giant overshoot" has been observed 
experimentally [8|, |2l|. Even an overload equal to the 
saturated stress deep in the silo produces a local rise in 
the stress instead of the flat stress profile predicted by the 
Janssen model. This can in part be explained by local 
changes in the mobilization of friction and packing frac- 
tion of granular matter induced by placing an overload 
on the material. 

To account for the overshoot effect and better fit ex- 
perimentally observed stress profiles, a piece-wise model 
has been proposed that describes stress in a silo as hy- 
drostatic - increasing linearly with depth - until a tun- 
able depth where the stress then becomes exponential 
(Tol [Tsj . This piece- wise behavior arises out of the anal- 
ysis for horizontally varying vertical stress and indicates 
the de pth to which the stress at the top boundary pene- 
trates [20[ . Fitting experimental results to the piece- wise 
pressure profile corrects experiments that find unphysi- 
cal values of the redirection parameter when the basic 
Janssen model is used 8, 10]. 

In some experiments, a threshold for overloads has 
been observed that, when passed, causes the redirection 
parameter k to increase [Ig] . When the presence of force 
networks within granular silos is considered, the mech- 
anism accounting for different overload profiles becomes 
clearer. Granular matter has been described as being in 
a "fragile" state. Along the directions in which there are 
pre-existing force chains, large loads can be supported. 
In other directions, only slight deformation of the inter- 
particle contacts are allowed before the material must 
rearrange, forming new force chains [2^ . To understand 
how a granular silo responds to overloads, we need to 
consider how the force network facilitates the redirection 
of stress to the walls, which is the goal of the present 
experiments. 

Since the Janssen model is a continuum description 
that at least moderately successfully predicts stress sat- 
uration with depth, it is interesting to examine in what 
ways it gives an appropriate mean description of the 
forces within a granular material. It is equally interesting 
to determine how and why it fails P, [2^. For instance, 
various models of granular force propagation give differ- 
ent predictions of mean silo behavior [2 [13, [1^ . Further 
experimental analysis of the silo system can serve to test 
these and future models of granular matter. 

In our experiments, we use photoelastic particles to 
investigate how the preparation and properties of the 
granular material and the network of contact forces in- 
teract to redirect stress to the walls and determine the 
development of the stress profile with depth. We con- 
duct experiments in a quasi-two-dimensional system that 
is described by the Janssen analysis above and has been 
shown to be similar to three dimensional systems in ex- 
periments and simulations ^3i!2ljll3l- Similar quasi- two- 
dimensional experiments have used photoelasticity in the 
past to examine the propagation of stress but did not test 
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FIG. 3: We examine photoelastic disks constrained to a quasi- 
two dimensional geometry in a vertical 'silo' made from two 
Plexiglas sheets separated by slightly more than a disk thick- 
ness. The bi-disperse disks are confined by vertical aluminum 
side walls to a channel L — 10 cm wide and rest upon a piston 
that can be used to control the mobilization of friction. The 
aspect ratio of the pile is typically 4:1. 

the Janssen model [1, [H, [2^ . Rather than measure the 
apparent weight at the bottom of a three-dimensional 
pile, these photoelastic techniques allow us to examine 
the actual force network responsible for Janssen screen- 
ing. We investigate the role of particle elasticity and force 
networks within silos by examining the mean response to 
overloads, by characterizing the properties of the force 
network, and by varying the elasticity of the particles. 



III. METHODOLOGY 

We study a quasi- two dimensional arrangement of pho- 
toelastic disks constrained vertically by two aluminum 
walls and on the faces by two transparent sheets to form 
a channel L = 10 cm wide, as depicted in Fig. [3) To pre- 
vent ordering, we use a bi-disperse mixture of disks, with 



the larger disks approximately 20% larger. The silo is 
constructed so that there are movable pistons at the top 
and bottom of the silo. The pile rests upon the bottom 
piston, while loads can be applied to the top piston. A 
lOg load applied to the top of the pile corresponds to a 
~ 1000 dyne/cm stress at the top. Both pistons can be 
pinned to prevent motion or unintended loading. 

We measure the local disk-scale pressure, which we 
refer to as "stress" on our photoelastic disks using a 
technique similar to previous granular research from our 
group 3]. The bi-refringent disks display bright fringes 
in response to applied stress when illuminated between 
crossed polarizers. The density of fringes within a given 
disk is proportional to the stress within the disk. We 
determine the stress on a given disk by first calculat- 
ing the gradient of an image of fringes. We take the 
square of the gradient averaged over multiple directions 
to be proportional to the number of fringes. We use 
an experimentally-obtained calibration to convert the 
squared gradient of the intensity into an equivalent load 
on a monolayer of disks and hence a single disk. 

We performed experiments with two different 
silo/camera arrangements. In the first silo we use 
approximately 650 disks — one quarter with diameter 
0.9 cm and the rest 0.75 cm — resting upon a piston 
that is affixed to a stepper motor. A 640 x 480 digital 
video camera is used to image the column at a rate of 30 
frames/s while the piston supporting the pile is slowly 
lowered at a rate of hundredths of a grain diameter per 
second. We examine the sequence of images to find 
the first frame where an abrupt rearrangement of the 
force network is observed and take the frame before that 
as the moment that the friction of the grains at the 
walls is maximally mobilized upwards. We perform our 
analysis on this maximally mobilized frame. To prepare 
for each observation, a new force network is created 
by quickly forcing the lower piston upwards and then 
letting the entire pile rain downwards as the piston falls. 
Simulations have shown that granular piles created by 
'raining' particles from above are independent of the rate 
at which grains are added [l^ . To apply an overload we 
gently lower the upper piston onto the top of the pile 
before the pile is lowered by the stepper motor. The 
top layers of grains rearrange slightly as the overload is 
applied but the remainder of the particles maintain their 
relative positions. 

In the second silo arrangement, approximately 1200 
slightly smaller disks are used (an approximately 3:1 mix- 
ture of 0.6 cm and 0.75 cm disks). The lower piston rests 
upon a base that can be manually lowered a small dis- 
tance. Without a fixed connection to a stepper motor, we 
are able to create new force networks by physically lifting 
the entire silo and hipping end-over-end. After the pile 
is lowered slightly to mobilize friction more uniformly, a 
high-resolution (3264 x 2468) digital still camera is used 
to take two pictures, one without polarization that can be 
used to find particle centers and one with polarization for 
obtaining force from photoelasticity. We can use this in- 
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FIG. 4: (Color Online) These composite images, generated by 
merging sets of seventy to one hundred gradient squared im- 
ages of force networks, display the onset of Janssen screening 
as the overload is increased from to 106 g in approximately 
25 g increments. The same intensity scale is used for all im- 
ages, indicating greater pressure deep in the silo for smaller 
loads. 

formation to characterize the contact and force networks 
within the material. To investigate the influence of parti- 
cle elasticity on these networks, we use disks made from 
two photoelastic polymers of different stiffness (Young's 
modulus E ^ A MPa for softer particles, E = 210 MPa 
for harder particles). Both soft and hard disks are 0.635 
cm thick, and the material from which they are made has 
a density of 1.2g/cm^. 

For all figures, the standard error indicated by bars 
for mean quantities is determined using a 'bootstrap' al- 
gorithm where the mean of the quantity is calculated 
for each of 200 randomized data sets constructed by se- 
lecting elements of the original data set at random with 
replacement [1^ . The standard deviation of these means 
is taken to be the standard error. 



IV. RESULTS 
A. Fixed Silo Method-Softer Particles 

We first used the fixed silo and stepper-motor arrange- 
ment to observe ^ 100 unloaded force networks and ~ 75 
networks each for four different overloads. Figure |4] shows 
the merged squared gradient images of these force net- 
works, roughly depicting mean response. If we bin the 
responses across the width of the silo into horizontal slices 
the height of one particle, we can obtain a pressure profile 
with depth as in Fig. [5l 

Considering first the case of the silo without an over- 
load, we see distinctly non-Janssen behavior. Although 
the pressure initially seems to saturate with moderate 
depth, the pressure starts increasing again toward the 



FIG. 5: By binning the images of force networks horizon- 
tally, averaging sets of seventy to one hundred images and 
then using a calibration to convert intensity gradient into 
pressure, we can generate pressure profiles. The dashed line 
indicates the unloaded case, while the increasing solid lines 
indicate overloads applied at the top increasing in increments 
of roughly 25 g; no load (•), 31 g (■), 56 g (a), 81 g (♦) 
and 106 g (T). Error bars indicate the standard deviation 
of means calculated for an ensemble of subsets of the data. 
Depth is given in terms of the width of the system L = 10 cm 
~ 13 smaller disk diameters. 



bottom of the silo. This increase is also observed for all 
overloads applied to the top of the pile. This effect is at 
least partly due to the way that the preparation history 
affects the unloaded state. As seen in the left-most im- 
age of Figure SI the mean stress is not constant across 
the layer at a given height, particularly near the bottom. 
This effect is presumably associated with the fact that 
when the piston pushes upwards, not all the particles are 
easily thrown upward, or that there is wall drag as the 
particles fall which affects their packing. 

We next examined the response to four overloads in in- 
crements of roughly 25 g (i.e. a two-dimensional stress of 
roughly 2500 dyne/cm). As shown in Fig.[5l for the two 
largest overloads a non-monotonic overshoot effect is ob- 
served. This behavior is similar to previous experiments 
[1, [2l[ - after an initial increase in pressure at the top 
of the silo, the largest overloads are screened by arching 
and the pressure converges with the unloaded pressure 
profile with depth. 

The behavior of the two smaller overloads is strikingly 
different from expected and previously-observed behav- 
ior. The additional pressure of small overloads persists 
with depth and is not screened by wall friction. The 
difference between large and small overloads is clear in 
Fig. [51 where the unloaded pressure profile has been sub- 
tracted from the four loaded pressure profiles. Deep 
within the pile, the additional pressure of two larger over- 
loads decays to zero, while the pressure from the smaller 
two loads persists. 
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FIG. 6: Plotting the difference between the pressure profile 
resulting from an overload and the unloaded profile shows 
how the overload is distributed. Four different overloads, 31 
g (■), 56 g (a), 81 g (♦) and 106 g (T), were applied across 
the top of the system. Smaller loads seem to penetrate deeper 
into the pile, while larger loads seem to show the previously 
observed overshoot effect followed by the expected decay. 
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FIG. 7: The Janssen length-scale A can be calculated nu- 
merically using the density of the granular material and the 
pressure profiles for no load (•), 56 g (a), 81 g (♦) and 106 g 
(T). We observe that A varies as a function of depth instead 
of remaining constant as in the Janssen model. Since the cal- 
culated length-scale varies with the inverse of the derivative 
of pressure, variations in the measured pressure profile can 
be greatly amplified, as with the for 31 g results at shallow 
depth. 



Given the non-monotonic, and sometimes non- 
saturating behavior of the observed pressure profiles, it 
is not surprising that the Janssen length-scale A — 2jt|— 
varies with depth, as indicated in Fig. [T] Using our data. 



FIG. 8: The Janssen screening length-scale A as a function 
of depth for, from top to bottom, no load (•), a 31 g load 
(A)and an 81 g (♦) load for both the initial (dashed) force 
network and the maximally mobilized (solid) force network. 



we estimate A numerically as: 



A(z) = 



P(z) 



Pps<l>g 



dP(z) 



(5) 



where is the density of the photoelastic material, (p is 
the packing fraction in the silo, P{z) is the horizontally- 
binned pressure and ^^^^ is estimated numerically us- 
ing local quadratic fits to P{z). Since A varies as the 
inverse of the redirection parameter fc, this result indi- 
cates that the assumption that k is constant with depth 
within our system is incorrect. Though the scatter in 
the pressure profile for 31 g introduces substantial noise 
into the calculation of A, hence this data is omitted, we 
see that the variation of A with depth is roughly similar 
for both unloaded and loaded silos. This indicates that 
the previously-observed load-induced variation in k is not 
responsible for variation in load propagation (l9| . 

Though Janssen's k varies with depth, we can still in- 
vestigate if k{z) behaves as a depth-dependent redirec- 
tion parameter might be expected to behave. In Fig. [H 
we calculate A as a function of depth for the ensemble of 
unloaded networks, and the ensembles for overloads just 
above and below the transition from deep-propagation 
to screening. We compare the network in the first im- 
age of our lowering sequence, corresponding to the initial 
network before the pile has been lowered to change mobi- 
lization, to the network in the subsequent image where we 
assume friction is maximally mobilized. Since A varies as 
the inverse of the mobilization k, the observed decrease in 
A after lowering the pile indicates that we have increased 
the mobilization. The mobilization does not change sig- 
nificantly at the top of the pile, indicating that our low- 
ering technique does not uniformly effect the pile. The 
difference between initial and mobilized frames is simi- 
lar for both small, deep-propagating overloads and large, 
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FIG. 9: By changing our procedure for creating new force 
networks we obtained pressure profiles that indicate pressure 
saturation with depth for an unloaded silo (•) and overloads 
of 20 g (■), 40 g (A), 60 g (♦) and 80 g (T). Despite the 
change in pile creation procedure the overshoot effect is still 
observed for large loads and small loads are still observed to 
penetrate to a considerable depth. 



decaying overloads. 

We note four important differences between the pres- 
sure profiles in Fig. [5] and those predicted by Janssen- 
like models, e.g. Fig. [5J A, and hence fc, is not constant; 
pressure from larger overloads first increases with depth 
before decaying; overall pressure does not seem to satu- 
rate with depth; and most strikingly, pressure from small 
overloads does not decay. The first two differences have 
been observed previously and are not surprising, given 
the previously discussed drawbacks of Janssen's analysis 
[1, [13, HH . To better understand the lack of overall pres- 
sure saturation, we conducted experiments in our second 
apparatus that allowed the particles to be rained from 
above, rather than pushed from the bottom. As shown 
in Fig. [5] we found that the second preparation procedure 
did create pressure profiles that saturate with depth. The 
largest overload is observed to overshoot before decaying 
and the smallest load propagates to great depth. 

Numerically determining the Janssen screening-length 
A for the saturating pressure profiles via Eq. O indicates 
that after an initial transient, the length-scale actually 
converges to a constant value of roughly 20 cm with 
depth, as in Fig. [TOl This implies a nominal value of 
= L/2A = 0.25. Interestingly, the \{z) for the moder- 
ate overloads are the flattest, while the screening length 
converges from a larger value for the larger overloads and 
from a smaller value for the silo with no overload and 
the smaller deep-penetrating 20 g overload. The non- 
constant value for the unloaded case may be due to sev- 
eral causes. One may be an underestimation of the small 
forces at the top of our pile when using our gradient 
method to calibrate photoelastic fringes. Alternatively, 
the variation at the top of the pile may be explained as 



FIG. 10: The numerically evaluated Janssen screening length- 
scale for the pressure profiles for an unloaded silo (•) and 
overloads of 20 g (■), 40 g (a), 60 g (♦) and 80 g (T) shows 
convergence with depth to a nearly constant value of A. 
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FIG. 11: For hard {E = 210 MPa) photoelastic disks we only 
see the overload profile due to the relative insensitivity of the 
harder material to stress. We examined overloads of 1837 g 
(•), 2267 g (■) and 3402 g (A). All seem well-described by 
exponential fits (dashed lines). 



a combination of both preparation — the act of placing 
an overload at the top — and non-uniform response to 
lowering the pile with depth. These two factors may also 
account for the giant overshoot effect seen for the large 
overload. 



B. Silo Flipping Method— Softer and Harder 
Particles 



We then used the silo flipping preparation to examine 
both the effect of preparation and of changing the elas- 
ticity of the particles, from Young's modulus E = 4 MPa 
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FIG. 12: Due to the insensitivity of the harder photoelas- 
tic particles, determining the local Janssen screening length- 
scale A numerically is only possible for the largest 3402 g (a) 
overload. The observed screening length is roughly half as 
long as for the soft {E — 4 MPa) particles. 



to E = 210 MPa. For the harder photoelastic material, 
the sensitivity to applied stress was lower, with the re- 
sult that the stress due to gravity was not visible even 
at great depth in the silo. In this case, any pressure pro- 
file is effectively an overload decay profile, since there is 
no unloaded profile to subtract. As shown in Fig. [TTl 
the additional pressure due to the overloads decays with 
depth after overshooting for all three cases. For all three 
overloads, the mean pressure eventually drops below the 
range of observable sensitivity, though large fluctuations 
about the mean were still present. All three cases can be 
fit with exponentials if the initial overshoot is neglected. 
Since the non-decaying, deep-penetrating overloads for 
the softer particles were close to the saturation pressure, 
it is not possible to determine whether or not such effects 
are present in the harder particles using our method. 

Determining the Janssen screening length-scale A nu- 
merically for the hard particle overloads was problematic 
since the pressure profiles eventually decayed below the 
sensitivity of the hard photoelastic material. In Fig. [T^ 
we show the calculated X{z) for the largest, 3402 g over- 
load in the top third of the silo. What we see is qual- 
itatively similar to the behavior of the largest overload 
in Fig. [TO] — A approaches a constant value from above. 
Surprisingly, although the system geometry is the same, 
the screening length is roughly half of what was found for 
the soft {E — A MPa) particles used to obtain the first 
two pressure profiles (Fig. O and Fig. [5]). 

The two sets of experiments nominally contrast soft 
and hard particles. However, there is another aspect 
which may be at least as important as granular elastic- 
ity. We note that the ratio of the applied loads to the 
Young's modulus, E, are comparable for the hard and 
the soft material. For the maximum load on the softer 
particles, ~ lOOg, this ratio is ~ 1.6 x 10~^ and for the 
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FIG. 13: Probability distributions for the difference in pres- 
sure from the mean value at the center (solid) and sides 
(open) of the silo for soft (top) particles in an unloaded silo 
(•) and with overloads of 20 g (■) and 80 g (a). We also show 
distributions for hard (bottom) particles with overloads of 
1837 g (•), 2267 g (■) and 3402 g (a). The curves differ for 
low pressures due to the differences in sensitivity of the hard 
and soft photoelastic polymers, but generally there is a great 
deal of overlap and we cannot distinguish the sides from the 
center. 



maximum load on the harder particles, 3400(7, this ra- 
tio is ^ 1.0 X 10^^. Thus, in appropriately scaled terms 
the loads are comparable in both cases, and hence on 
this basis, elasticity cannot explain the difference in the 
experimental observations. However, the ratio of the ap- 
plied load to the weight of the material is significantly 
different in the two cases, namely about 0.5 for the soft 
particles, and 17 for the hard particles. Thus, we con- 
clude that for the softer particles, the loads are not large 
enough to overcome the effects caused by gravitational 
loading, which necessarily reflects the preparation his- 
tory. 
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C. Horizontal Variations and Fluctuations 

It is also interesting to investigate various other aspects 
of the force profiles which differ from the usual Janssen 
picture. In particular, we consider horizontal variabil- 
ity of the forces, and also force fluctuations. In this re- 
gard, we note model studies by Pitman which exhibited 
horizontal de pen dence on fluctuations and stress within 
granular silos [l6| due to inhomogeneous mobilization of 
friction. 

As shown in Fig. [131 we measured the local deviation 
from the mean pressure at a given height for both the 
center and sides of the silo for all overloads for both hard 
and soft particles. We did not see significant horizontal 
variation in the distributions, and for the soft particles 
the distributions were nearly identical for all cases and 
very reminiscent of force distributions observed for other 
granular systems [13] ■ Forces below the mean fall off 
quickly, while forces larger than the mean tail off slowly, 
and roughly follow an exponential. 

Recent simulations for vertical granular matter with 
periodic boundary conditions (in the horizontal direc- 
tions) have indicated that the averaging length-scale suf- 
ficient to determine the mean pressure is possibly as small 
as a grain diameter [3l| . Though our geometry has side 
walls, it is still worthwhile to determine the length over 
which we must average for the pressure profile to con- 
verge. To do this, we determine a grid of points within 
the silo and average the pressure within a tunable radius 
to assign a pressure to that point. In Fig.[T31 we plot the 
pressure in the middle of the silo (on axis) with depth as 
determined by averaging over regions with radii of 1, 5, 
10 and 15 small particle diameters. For the two larger 
averaging regions, the profiles converge to those observed 
using horizontal binning by depth in Fig. [9l This con- 
vergence indicates a correlation length-scale between 10 
and 15 particle diameters — approximately the same as 
the Janssen screening length-scale obtained numerically. 
This is interesting since in this instance, the macroscopic 
scale of the Janssen parameter A is comparable to the 
averaging scale, even though the Janssen model does not 
include particle size. As can be seen for the lower image 
in Fig. [TU the the effect of averaging length is roughly 
the same regardless of particle stiffness. 

An additional aspect of interest is the geometric nature 
of the contacts. Relevant quantities include the packing 
fraction, 0, the number of contacts per particle, Z, and 
clustering coefficient, defined below. 

With our second preparation procedure, we are able 
to obtain to obtain all of these quantities. We do this 
by capturing both polarized images showing the force 
network within the photoelastic disks and unpolarized 
images showing the locations of those particles. Particle 
centers can be identified using a convolution with sep- 
arate kernels for the two particle sizes. We then use a 
numerical algorithm to reject overlapping particles and 
identify neighbor particles in contact with each particle 
based upon the geometric consideration that the distance 
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FIG. 14: Pressure profiles on axis converge as tfie length R 
over whicfi averaging is conducted varies, as indicated by lines 
of increasing thickness for R = 1, 5, 10 and 15 small particle 
diameters for both soft particles (top) in an unloaded silo (•) 
and with overloads of 20 g (■) and 80 g (a) and for hard 
particles (bottom) with overloads of 1837 g (•), 2267 g (■) 
and 3402 g (a). 



between two particles must be equal to the sum of their 
radii. This allows us to determine a local packing fraction 
throughout the silo. 

Since previous studies have shown that the packing 
fraction of a pile relates to the redirection parameter, 
we investigate depth-dependent variations of the packing 
fraction [8| . Due to the dilation of the material, we find 
that the packing fraction is slightly lower for loaded silos 
for the depths available in the unloaded case. We do not 
find systematic correlation between the packing fraction 
and the screening length-scale. 

The particle locations and contacts comprise a contact 
network for each pile. Analyzing these contact networks 
allows the determination of graph and network properties 
such as the average number of contacts per particle, Z. 
We can also determine the 'clustering coefficient', which 
gives a measure of how tightly interconnected is the net- 
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FIG. 15: The average packing fraction for soft (left) and 
hard (right) particles varies slightly depending upon over- 
load. While it is roughly constant with depth for an unloaded 
silo (•) it varies for overloads of 20 g (■) and 80 g (a). 



work 32]. For each particle we consider the n neighbor- 
ing particles with which that particle is in contact. There 
are n — 1 possible neighbor-neighbor contacts for round 
particles in two dimensions, so we define the clustering 
coefficient as the ratio of the m actual contacts present 
between those neighbors to the n — 1 possible contacts, 

71—1 

For each contact network, we can determine a subset 
of contacts for which the particle-scale pressure is greater 
than the mean for the whole silo. This subset of contacts 
can be thought of as roughly comprising a force network, 
for which the same graph and network pro perties as the 
contact network can also be calculated [33| . 

We first examine the dependence of the mean number 
of contacts per particle with depth for both the over- 
all contact network and the force network. As shown in 
Fig. [161 we find that for the contact network, the number 
of contacts is nearly independent of overload and con- 
stant with depth at a value close to A; = 4, the typically 
observed coordination number for a disordered granular 
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FIG. 16; The average number of contacts per particle is 
roughly constant for the contact network (left) but varies 
at the top of the pile for the force network (right). This is 
true for both soft particles (above) in an unloaded silo (•) 
and with overloads of 20 g (■) and 80 g (k) as a well as for 
hard particles (below) with overloads of 1837 g (•), 2267 g 
(■) and 3402 g (a). 



material, although slightly lower due to the presence of 
the walls [S^. The average number of contacts in the 
force network is slightly less than to two, which is not 
surprising since many particles in a chain have only two 
contacts and, in our relative small silo, many particles 
make contact with the walls. 

Significantly, for the force network within both the 
hard and soft particles, the average number of contacts 
at the top of the pile varies depending upon the over- 
load, before converging with depth, unlike the contact 
network which does not vary with overload. This indi- 
cates that the force network is quantitatively different, 
depending upon the size of the applied overload. For the 
deep-penetrating overload on the soft particles, the mean 
number of contacts at the top is closest to the saturation 
value, perhaps indicating that the structure of the force 
network is made visible by the additional pressure, but 
is not disturbed as in the case of larger overloads. 

For the clustering coefficient, Fig.fTTlshows the contact 
network to be better connected than the force network. 
As with the number of contacts, the clustering coefficient 
is roughly constant for the contact network but varies 
with depth for the force network. For hard particles, the 
overloads do not seem to alter the clustering coefficient 
of the contact or force network. 

The differences between the contact and force networks 
may indicate an explanation for the deep propagation 
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FIG. 17: For the contact network (left) in soft (above) and 
hard (below) particles the average clustering coefficient is 
roughly the same for an unloaded silo (•) and overloads of 
20 g (■) and 80 g (A) for soft particles and 1837 g (•), 2267 
g (■) and 3402 g (a) for hard particles. For the force net- 
work (right) within soft particles, the clustering coefficient 
varies depending upon overload. For hard particles the force 
network is unchanged by overload. 



of small overloads. If the material within a silo can be 
characterized as fragile matter, then we expect that it will 
have the ability to support forces along the force chains 
that are present, but that the force network will have 
to rearrange to support forces along normal directions 
[2^ or to carry large loads. This is borne out for our 
softer particles which are more photoelastically sensitive. 
We observe deep-penetrating response for small loads for 
which there is only slight deformation of the network, 
while larger loads rearrange the force network, resulting 
in Janssen screening. 



V. CONCLUSION 

For elastic particles in granular silos, we observe a non- 
linear response as we vary the size of overloads placed 
upon the material at the top of the silo. For sufficiently 
small overloads, the pre-existing contact network con- 



veys additional pressure deep into the pile, contrary to 
the classical Janssen description. For larger overloads we 
find that the force network rearranges, creating frictional 
arches braced against the walls that act to screen the 
depths of the silo from the additional pressure. Though 
the Janssen model predicts monotonic saturation, we find 
that the pressure initially increases with depth before de- 
caying — the giant overshoot phenomenon ^] . 

By numerically evaluating the Janssen length-scale we 
find that it is not constant, meaning that the simple in- 
tegration typically used to derive pressure profile predic- 
tions from force balance is inadequate in our experiments, 
and therefore quite likely in general. Wc find that there 
is a different Janssen length-scale at the top of the pile, 
possibly corresponding to changes in mobilization result- 
ing from placing the overload onto the material. 

We find that the pressure profile is sensitive to prepa- 
ration, and that this sensitivity depends on both loading 
and particle stiffness. The greatest dependence on his- 
tory occurs when the magnitudes of gravitational loading 
and the applied loading are commensurate. This situa- 
tion occurs for our softer particles, for which we apply 
a maximum load that corresponds to roughly half the 
weight of the particles. For our harder particles, which 
have weight comparable to our softer particles, all ap- 
plied loads significantly exceed the particle weight. Such 
loads lead to particle deformations comparable to the soft 
case, but the observed non-linear behavior of the profiles 
is substantially less. 

We find that distributions for the particle-scale pres- 
sure are similar to those observed for local force measures 
in other granular systems. Specifically, there is a roughly 
exponential decay of the distributions at large pressures. 
There is a peak in the distributions near the mean, and 
a rapid fall off below the mean. In general, we find that 
force networks exhibit fluctuations with a length-scale 
comparable to the width of the silo. 
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